"lost.job12:cvd.gmb.gvt" = "Kept Job, CV-19\n(Baseline = Lost Job):Trust Gov., CV-19",
"lost.job1:cvd.gmb.gvt" = "Lost Job, CV-19\n(Baseline = No Job):Trust Gov., CV-19",
"lost.job2:cvd.gmb.gvt" = "Kept Job, CV-19\n\n(Baseline = No Job):Trust Gov., CV-19")
## Main Analyses: Aspiration -- Including Covid Daily Cases and Clustered SEs -- with interactions
texreg(list(sel.mg$b, sel.mg$beff, sel.mg$beffg,sel.mg$bgmigdems, sel.mg$bgmigdemsint0), custom.coef.map = keep.vars,
booktabs = T, dcolumn = T, stars = c(0.01, 0.05, 0.1), caption = "Main Analyses: Aspiration to Migrate", caption.above = T,
custom.note = "SEs clustered by region")
int.blnojob.factor <- lm(mg.asp ~ factor(lost.job) + cvd.hhs + cvd.cases.day
+ cvd.gmb.vacc.rv + cvd.gmb.econ.rv + cvd.rch.vacc.rv + cvd.rch.econ.rv
+ cvd.gmb.gvt + mg.ab + rskavrs + nt.mg.out
+ male + age + educ + factor(lost.job)*cvd.gmb.gvt, data=bl.f)
##generating clustered standard errors
int.blnojob.factor.clust <- cluster.vcov(int.blnojob.factor, bl.f$region.cln)
int.blnojob.factor.mcse <- coeftest(int.blnojob.factor, int.blnojob.factor.clust)
int.blnojob.factor.mcse
## Pull out data for marginal effects plots
# coefficients
int.blnojob.factor.b2 <- int.blnojob.factor.mcse[2,1]
int.blnojob.factor.b3 <- int.blnojob.factor.mcse[10,1]
int.blnojob.factor.b22 <- int.blnojob.factor.mcse[17,1] # interaction term
# variances
int.blnojob.factor.varb2 <- int.blnojob.factor.clust[2,2]
int.blnojob.factor.varb3 <- int.blnojob.factor.clust[10,10]
int.blnojob.factor.varb22 <- int.blnojob.factor.clust[17,17]
# covariances
int.blnojob.factor.covb2b22 <- int.blnojob.factor.clust[2,17]
int.blnojob.factor.covb3b22 <- int.blnojob.factor.clust[10,17]
##transforming dataframe to graph
int.blnojob.factor.df <- bl.f %>% mutate(cvd.gmb.gvt= seq(1, 5, length.out=829),
int.blnojob.factor.conbx = int.blnojob.factor.b2+int.blnojob.factor.b22*cvd.gmb.gvt,
int.blnojob.factor.consx = sqrt(int.blnojob.factor.varb2 + int.blnojob.factor.varb22*(cvd.gmb.gvt^2)+
2*int.blnojob.factor.covb2b22*cvd.gmb.gvt),
int.blnojob.factor.ax = 1.96*int.blnojob.factor.consx,
int.blnojob.factor.upperx = int.blnojob.factor.conbx + int.blnojob.factor.ax,
int.blnojob.factor.lowerx = int.blnojob.factor.conbx - int.blnojob.factor.ax)
# a rug
gmb.gvt.temp <- bl.f$cvd.gmb.gvt
a<-runif(829,min=0,max=0.6)
a<-a-0.2
gmb.gvt.jitter <-gmb.gvt.temp+a
blnojob.factor.plot <- ggplot(data = int.blnojob.factor.df) +
scale_y_continuous("Marginal Effect of Losing Job\n(Baseline: No Job)") +
scale_x_continuous(limits = c(1,5)) +
geom_ribbon(aes(x=cvd.gmb.gvt, ymin=int.blnojob.factor.lowerx, ymax = int.blnojob.factor.upperx), alpha=.2) +
geom_line(aes(x=cvd.gmb.gvt, y=int.blnojob.factor.conbx, color="Lost Job"), size = 1)  +
geom_hline(yintercept = 0) +
geom_rug(aes(x=gmb.gvt.jitter),sides="b") +
scale_fill_manual(values = wes_palette(n=2, name="Darjeeling1"),
name = "",
labels = c("Lost Job")) +
xlab("Trust Gambian Government to Handle Pandemic") +
ggtitle("Marginal Effect of Losing One's Job\n on the Aspiration to Migrate")+
theme_bw(base_size = 12)+
theme(text = element_text(size=20),
legend.position="none",
plot.title = element_text(face="bold", hjust = 0.5))
blnojob.factor.plot
texreg(list(sel.mg$bgmigdemsj0,sel.mg$bgmigdemsj1,sel.mg$bgmigdemsj2), custom.coef.map = keep.vars,
booktabs = T, dcolumn = T, stars = c(0.01, 0.05, 0.1), caption = "Different Baselines: Aspiration to Migrate", caption.above = T)
f.ol <- as.formula(paste("mg.asp.f", paste(mg.vars.f, collapse = " + "), sep = " ~ "))
b.ol <- polr(f.ol, data = bl.f, Hess = T)
b.ol.cl <- coeftest(b.ol, vcov=vcovCL(b.ol, factor(bl.f$region.cln) ))
texreg(list(sel.mg$bgmigdems, b.ol.cl), custom.coef.map = keep.vars,
booktabs = T, dcolumn = T, stars = c(0.01, 0.05, 0.1), caption = "Ord Log: Aspiration to Migrate", caption.above = T,
custom.model.names = c("OLS","Ord.L"))
cvd.att <- c("cvd.gmb.vacc", "cvd.rch.vacc","cvd.gmb.econ", "cvd.rch.econ")
cvd.df <- bl.f[,c("OLDUID",cvd.att)]
long.cvd <- melt(cvd.df, id.vars = c("OLDUID"), measure.vars = cvd.att, value.name = "timelines", na.rm = T)
long.cvd.vacc <- subset(long.cvd, long.cvd$variable == "cvd.gmb.vacc" | long.cvd$variable == "cvd.rch.vacc" )
long.cvd.econ <- subset(long.cvd, long.cvd$variable == "cvd.gmb.econ" | long.cvd$variable == "cvd.rch.econ" )
pd <- position_dodge(0.2) # move them .05 to the left and right
p1 <- ggplot(data=long.cvd.vacc, aes(x=timelines, group=variable, color=variable, fill = variable)) +
geom_density(adjust=2, alpha=0.02, na.rm = T, size = 0.7) +
ylab("") + labs(title="Vaccine Distribution") +
scale_color_manual(values = wes_palette(n=2, name="Darjeeling1"),
name = "",
labels = c("The Gambia", "Rich Countries")) +
scale_fill_manual(values = wes_palette(n=2, name="Darjeeling1"),
name = "",
labels = c("The Gambia", "Rich Countries")) +
scale_x_discrete("", limits = c("Already", "Within year","Within 2 Years", "Distant future", "Never"), expand = c(0.03,0.03)) +
scale_y_continuous(limits = c(0,1), expand = c(0,0)) +
theme_bw(base_size = 17) + theme(plot.title = element_text(size = 17))
p1
p2 <- ggplot(data=long.cvd.econ, aes(x=timelines, group=variable, color=variable, fill = variable)) +
geom_density(adjust=2, alpha=0.02, na.rm = T, size = 0.7) +
ylab("") + labs(title="Economic Recovery") +
scale_color_manual(values = wes_palette(n=2, name="Darjeeling1"),
name = "",
labels = c("The Gambia", "Rich Countries")) +
scale_fill_manual(values = wes_palette(n=2, name="Darjeeling1"),
name = "",
labels = c("The Gambia", "Rich Countries")) +
scale_x_discrete("", limits = c("Already", "Within year","Within 2 Years", "Distant future", "Never"), expand = c(0.03,0.03)) +
scale_y_continuous(limits = c(0,1.1), expand = c(0,0)) +
theme_bw(base_size = 17) + theme(plot.title = element_text(size = 17))
p2
lj <- table(bl.f$lost.job, useNA = "ifany")
male <- table(bl.f$male, useNA = "ifany")
reg <- table(bl.f$region.cln, useNA = "ifany")
# Others
describe(bl.f$cvd.hhs)
describe(bl.f$cvd.gmb.gvt)
describe(bl.f$rskavrs)
describe(bl.f$nt.mg.out)
describe(bl.f$age)
describe(bl.f$mg.asp)
describe(bl.f$mg.ab)
describe(bl.f$educ)
describe(bl.f$cvd.gmb.vacc.rv)
describe(bl.f$cvd.gmb.econ.rv)
describe(bl.f$cvd.rch.vacc.rv)
describe(bl.f$cvd.rch.econ.rv)
prop.table(table(bl.f$cvd.gmb.vacc, useNA = "ifany"))
prop.table(table(bl.f$cvd.gmb.econ, useNA = "ifany"))
prop.table(table(bl.f$cvd.rch.econ, useNA = "ifany"))
lj <- table(bl.f$lost.job, useNA = "ifany")
male <- table(bl.f$male, useNA = "ifany")
reg <- table(bl.f$region.cln, useNA = "ifany")
# Others
describe(bl.f$cvd.hhs)
describe(bl.f$cvd.gmb.gvt)
describe(bl.f$rskavrs)
describe(bl.f$nt.mg.out)
describe(bl.f$age)
describe(bl.f$mg.asp)
describe(bl.f$mg.ab)
describe(bl.f$educ)
describe(bl.f$cvd.gmb.vacc.rv)
describe(bl.f$cvd.gmb.econ.rv)
describe(bl.f$cvd.rch.vacc.rv)
describe(bl.f$cvd.rch.econ.rv)
prop.table(table(bl.f$cvd.gmb.vacc, useNA = "ifany"))
prop.table(table(bl.f$cvd.gmb.econ, useNA = "ifany"))
prop.table(table(bl.f$cvd.rch.econ, useNA = "ifany"))
mg.vars.f
bl.covid.df <- bl.f[mg.vars.f]
names(bl.covid.df)
dim(bl.f)
dim(bl.covid.df)
bl.covid.df$lost.job <- as.numeric(bl.covid.df$lost.job)
corr <- round(cor(bl.covid.df, use = "na.or.complete"), 2)
corr
library(coefplot)
#install.packages("dotwhisker")
library(dotwhisker)
#install.packages("broom")
library(broom)
library(dplyr)
#install.packages("tidyr")
require(tidyr)
#install.packages("tidyverse")
require(tidyverse)
require(ggplot2)
require(stargazer)
#install.packages("huxtable")
require(huxtable)
library(Hmisc)
#install.packages("vctrs")
library(vctrs)
library(reshape2)
#install.packages("Rmisc")
library(Rmisc)
library(RColorBrewer)
#install.packages("wesanderson")
library(wesanderson)
#install.packages("doBy")
library(doBy)
#install.packages("binom")
library(binom)
#install.packages("miceadds")
library(miceadds)
library(gridExtra)
#install.packages("ggpubr")
library(ggpubr)
library(survey)
library(psych)
library(estimatr)
require(MASS)
require(Hmisc)
require(reshape2)
library(lmtest)
library("sandwich")
library(texreg)
library(interplot)
#devtools::install_github("strengejacke/strengejacke")
library(sjPlot)
#install.packages("ggeffects")
library(ggeffects)
#install.packages("ggcorrplot")
library("ggcorrplot")
#install.packages("Cairo")
library(Cairo)
setwd("~/Dropbox/DFID_IOM/Analysis/II_Covid_Replication")
load(file = "bl.f.RData")
mg.vars.b <- c("lost.job","cvd.hhs","cvd.cases.day")
mg.vars.e <- c(mg.vars.b,"cvd.gmb.vacc.rv","cvd.gmb.econ.rv","cvd.rch.vacc.rv","cvd.rch.econ.rv")
mg.vars.t <- c(mg.vars.e,"cvd.gmb.gvt")
mg.vars.m <- c(mg.vars.t,"mg.ab","rskavrs","nt.mg.out")
mg.vars.f <- c(mg.vars.m,"male","age","educ")
## Different baselines for job (Appendix C.1.)
full.nj <- mg.vars.f[mg.vars.f != "lost.job"]
mg.vars.bl0 <- c("lost.job",full.nj)
mg.vars.bl1 <- c("lost.job1",full.nj)
mg.vars.bl2 <- c("lost.job2",full.nj)
## Interaction between jobs and confidence in government (Main analyses)
mg.vars.int0 <- c("lost.job*cvd.gmb.gvt","lost.job",full.nj)
## List of all variable lists for OLS
mlist.ols <- list(mg.vars.b, mg.vars.e,mg.vars.t,mg.vars.m,mg.vars.f,mg.vars.bl0,mg.vars.bl1,mg.vars.bl2,mg.vars.int0)
## Object for storing regression results
sel.mg <- vector(length=length(mlist.ols), mode = "list")
names(sel.mg) <- c("b","beff","beffg","bgmig","bgmigdems","bgmigdemsj0","bgmigdemsj1",
"bgmigdemsj2","bgmigdemsint0")
## Running all OLS models
for (i in 1: length(mlist.ols))
{
f <- as.formula(paste("mg.asp",
paste(mlist.ols[[i]], collapse = " + "),
sep = " ~ "))
b.y.cl <- lm.cluster(f, cluster = "region.cln", data = bl.f) # Cluster standard errors by region
b.y.x <- texreg::extract(b.y.cl) # For using with texreg
sel.mg[[i]] <- b.y.x
}
## Variables and labels for regression table
keep.vars <- list("lost.job1" = "Lost Job, CV-19\n(Baseline = No Job)", "lost.job2" = "Kept Job, CV-19\n(Baseline = No Job)",
"lost.job10" = "No Job, CV-19\n(Baseline = Lost Job)", "lost.job12" = "Kept Job, CV-19\n(Baseline = Lost Job)",
"lost.job20" = "No Job, CV-19\n(Baseline = Kept Job)", "lost.job21" = "Lost Job, CV-19\n(Baseline = Kept Job)",
"cvd.hhs" = "HH Savings, CV-19", "cvd.gmb.vacc.rv" = "Effcy. CV-19 Vaccine, Gambia",
"cvd.gmb.econ.rv" = "Effcy. CV-19 Econ, Gambia", "cvd.rch.vacc.rv" = "Effcy. CV-19 Vaccine, Rich",
"cvd.rch.econ.rv" = "Effcy. CV-19 Econ, Rich", "cvd.gmb.gvt" = "Trust Gov., CV-19",
"rskavrs" = "Risk Averse", "mg.ab" = "Ability to Migrate",
"nt.mg.out" = "Migrant Networks",
"male" = "Male", "age" = "Age", "educ" = "Education","cvd.cases.day" = "CV-19 Cases Cum. PMP",
"lost.job20:cvd.gmb.gvt" = "No Job, CV-19\n(Baseline = Kept Job):Trust Gov., CV-19",
"lost.job21:cvd.gmb.gvt" = "Lost Job, CV-19\n(Baseline = Kept Job):Trust Gov., CV-19",
"lost.job10:cvd.gmb.gvt" = "No Job, CV-19\n(Baseline = Lost Job):Trust Gov., CV-19",
"lost.job12:cvd.gmb.gvt" = "Kept Job, CV-19\n(Baseline = Lost Job):Trust Gov., CV-19",
"lost.job1:cvd.gmb.gvt" = "Lost Job, CV-19\n(Baseline = No Job):Trust Gov., CV-19",
"lost.job2:cvd.gmb.gvt" = "Kept Job, CV-19\n\n(Baseline = No Job):Trust Gov., CV-19")
## Main Analyses: Aspiration -- Including Covid Daily Cases and Clustered SEs -- with interactions
texreg(list(sel.mg$b, sel.mg$beff, sel.mg$beffg,sel.mg$bgmigdems, sel.mg$bgmigdemsint0), custom.coef.map = keep.vars,
booktabs = T, dcolumn = T, stars = c(0.01, 0.05, 0.1), caption = "Main Analyses: Aspiration to Migrate", caption.above = T,
custom.note = "SEs clustered by region")
# running the model
int.blnojob.factor <- lm(mg.asp ~ factor(lost.job) + cvd.hhs + cvd.cases.day
+ cvd.gmb.vacc.rv + cvd.gmb.econ.rv + cvd.rch.vacc.rv + cvd.rch.econ.rv
+ cvd.gmb.gvt + mg.ab + rskavrs + nt.mg.out
+ male + age + educ + factor(lost.job)*cvd.gmb.gvt, data=bl.f)
##generating clustered standard errors
int.blnojob.factor.clust <- cluster.vcov(int.blnojob.factor, bl.f$region.cln)
int.blnojob.factor.mcse <- coeftest(int.blnojob.factor, int.blnojob.factor.clust)
int.blnojob.factor.mcse
## Pull out data for marginal effects plots
# coefficients
int.blnojob.factor.b2 <- int.blnojob.factor.mcse[2,1]
int.blnojob.factor.b3 <- int.blnojob.factor.mcse[10,1]
int.blnojob.factor.b22 <- int.blnojob.factor.mcse[17,1] # interaction term
# variances
int.blnojob.factor.varb2 <- int.blnojob.factor.clust[2,2]
int.blnojob.factor.varb3 <- int.blnojob.factor.clust[10,10]
int.blnojob.factor.varb22 <- int.blnojob.factor.clust[17,17]
# covariances
int.blnojob.factor.covb2b22 <- int.blnojob.factor.clust[2,17]
int.blnojob.factor.covb3b22 <- int.blnojob.factor.clust[10,17]
##transforming dataframe to graph
int.blnojob.factor.df <- bl.f %>% mutate(cvd.gmb.gvt= seq(1, 5, length.out=829),
int.blnojob.factor.conbx = int.blnojob.factor.b2+int.blnojob.factor.b22*cvd.gmb.gvt,
int.blnojob.factor.consx = sqrt(int.blnojob.factor.varb2 + int.blnojob.factor.varb22*(cvd.gmb.gvt^2)+
2*int.blnojob.factor.covb2b22*cvd.gmb.gvt),
int.blnojob.factor.ax = 1.96*int.blnojob.factor.consx,
int.blnojob.factor.upperx = int.blnojob.factor.conbx + int.blnojob.factor.ax,
int.blnojob.factor.lowerx = int.blnojob.factor.conbx - int.blnojob.factor.ax)
# a rug
gmb.gvt.temp <- bl.f$cvd.gmb.gvt
a<-runif(829,min=0,max=0.6)
a<-a-0.2
gmb.gvt.jitter <-gmb.gvt.temp+a
blnojob.factor.plot <- ggplot(data = int.blnojob.factor.df) +
scale_y_continuous("Marginal Effect of Losing Job\n(Baseline: No Job)") +
scale_x_continuous(limits = c(1,5)) +
geom_ribbon(aes(x=cvd.gmb.gvt, ymin=int.blnojob.factor.lowerx, ymax = int.blnojob.factor.upperx), alpha=.2) +
geom_line(aes(x=cvd.gmb.gvt, y=int.blnojob.factor.conbx, color="Lost Job"), size = 1)  +
geom_hline(yintercept = 0) +
geom_rug(aes(x=gmb.gvt.jitter),sides="b") +
scale_fill_manual(values = wes_palette(n=2, name="Darjeeling1"),
name = "",
labels = c("Lost Job")) +
xlab("Trust Gambian Government to Handle Pandemic") +
ggtitle("Marginal Effect of Losing One's Job\n on the Aspiration to Migrate")+
theme_bw(base_size = 12)+
theme(text = element_text(size=20),
legend.position="none",
plot.title = element_text(face="bold", hjust = 0.5))
blnojob.factor.plot
texreg(list(sel.mg$bgmigdemsj0,sel.mg$bgmigdemsj1,sel.mg$bgmigdemsj2), custom.coef.map = keep.vars,
booktabs = T, dcolumn = T, stars = c(0.01, 0.05, 0.1), caption = "Different Baselines: Aspiration to Migrate", caption.above = T)
f.ol <- as.formula(paste("mg.asp.f", paste(mg.vars.f, collapse = " + "), sep = " ~ "))
b.ol <- polr(f.ol, data = bl.f, Hess = T)
b.ol.cl <- coeftest(b.ol, vcov=vcovCL(b.ol, factor(bl.f$region.cln) ))
texreg(list(sel.mg$bgmigdems, b.ol.cl), custom.coef.map = keep.vars,
booktabs = T, dcolumn = T, stars = c(0.01, 0.05, 0.1), caption = "Ord Log: Aspiration to Migrate", caption.above = T,
custom.model.names = c("OLS","Ord.L"))
### Timeline Density Plots (Figure 1, main text)
cvd.att <- c("cvd.gmb.vacc", "cvd.rch.vacc","cvd.gmb.econ", "cvd.rch.econ")
cvd.df <- bl.f[,c("OLDUID",cvd.att)]
long.cvd <- melt(cvd.df, id.vars = c("OLDUID"), measure.vars = cvd.att, value.name = "timelines", na.rm = T)
long.cvd.vacc <- subset(long.cvd, long.cvd$variable == "cvd.gmb.vacc" | long.cvd$variable == "cvd.rch.vacc" )
long.cvd.econ <- subset(long.cvd, long.cvd$variable == "cvd.gmb.econ" | long.cvd$variable == "cvd.rch.econ" )
pd <- position_dodge(0.2) # move them .05 to the left and right
p1 <- ggplot(data=long.cvd.vacc, aes(x=timelines, group=variable, color=variable, fill = variable)) +
geom_density(adjust=2, alpha=0.02, na.rm = T, size = 0.7) +
ylab("") + labs(title="Vaccine Distribution") +
scale_color_manual(values = wes_palette(n=2, name="Darjeeling1"),
name = "",
labels = c("The Gambia", "Rich Countries")) +
scale_fill_manual(values = wes_palette(n=2, name="Darjeeling1"),
name = "",
labels = c("The Gambia", "Rich Countries")) +
scale_x_discrete("", limits = c("Already", "Within year","Within 2 Years", "Distant future", "Never"), expand = c(0.03,0.03)) +
scale_y_continuous(limits = c(0,1), expand = c(0,0)) +
theme_bw(base_size = 17) + theme(plot.title = element_text(size = 17))
p1
p2 <- ggplot(data=long.cvd.econ, aes(x=timelines, group=variable, color=variable, fill = variable)) +
geom_density(adjust=2, alpha=0.02, na.rm = T, size = 0.7) +
ylab("") + labs(title="Economic Recovery") +
scale_color_manual(values = wes_palette(n=2, name="Darjeeling1"),
name = "",
labels = c("The Gambia", "Rich Countries")) +
scale_fill_manual(values = wes_palette(n=2, name="Darjeeling1"),
name = "",
labels = c("The Gambia", "Rich Countries")) +
scale_x_discrete("", limits = c("Already", "Within year","Within 2 Years", "Distant future", "Never"), expand = c(0.03,0.03)) +
scale_y_continuous(limits = c(0,1.1), expand = c(0,0)) +
theme_bw(base_size = 17) + theme(plot.title = element_text(size = 17))
p2
# Proportions
lj <- table(bl.f$lost.job, useNA = "ifany")
male <- table(bl.f$male, useNA = "ifany")
reg <- table(bl.f$region.cln, useNA = "ifany")
# Others
describe(bl.f$cvd.hhs)
describe(bl.f$cvd.gmb.gvt)
describe(bl.f$rskavrs)
describe(bl.f$nt.mg.out)
describe(bl.f$age)
describe(bl.f$mg.asp)
describe(bl.f$mg.ab)
describe(bl.f$educ)
describe(bl.f$cvd.gmb.vacc.rv)
describe(bl.f$cvd.gmb.econ.rv)
describe(bl.f$cvd.rch.vacc.rv)
describe(bl.f$cvd.rch.econ.rv)
prop.table(table(bl.f$cvd.gmb.vacc, useNA = "ifany"))
prop.table(table(bl.f$cvd.gmb.econ, useNA = "ifany"))
prop.table(table(bl.f$cvd.rch.econ, useNA = "ifany"))
bl.covid.df <- bl.f[mg.vars.f]
dim(bl.f)
dim(bl.covid.df)
bl.covid.df$lost.job <- as.numeric(bl.covid.df$lost.job)
corr <- round(cor(bl.covid.df, use = "na.or.complete"), 2)
corr
colnames(corr) <- c("Lost Job", "HH Savings", "CV-19 Cases Cum", "Effcy. CV-19 Vaccine, Gambia",
"Effcy. CV-19 Econ, Gambia", "Effcy. CV-19 Vaccine, Rich", "Effcy. CV-19 Econ, Rich",
"Trust Gov., CV-19", "Ability to Migrate", "Risk Averse", "Migrant Networks",
"Male", "Age", "Education")
rownames(corr) <- c("Lost Job", "HH Savings", "CV-19 Cases Cum", "Effcy. CV-19 Vaccine, Gambia",
"Effcy. CV-19 Econ, Gambia", "Effcy. CV-19 Vaccine, Rich", "Effcy. CV-19 Econ, Rich",
"Trust Gov., CV-19", "Ability to Migrate", "Risk Averse", "Migrant Networks",
"Male", "Age", "Education")
ggcorrplot(corr,lab = TRUE, lab_size = 2, tl.cex = 7)
ggsave("corrmatrix.eps", width=12, height=10, unit="cm", device = cairo_ps)
cvd.cases <- read.csv("owid-covid-data_gmb_jan2021.csv")
cvd.cases$date2 <- as.Date(parse_date_time(cvd.cases$date,"mdy"))
cvd1 <- cvd.cases[which(cvd.cases$date2 > as.Date("2021-01-06") & cvd.cases$date2 < as.Date("2021-01-28")),]
png(filename = "Cvd_Cases.png", height = 480, width = 480)
p<-ggplot(data=cvd1, aes(x=date2, y=new_cases)) +
geom_bar(stat="identity") + labs(x= "Day", y="New Cases, Raw") +
scale_x_date(date_breaks = "1 day") + theme_bw(base_size = 12) +
theme(axis.text.x = element_text(angle = 90))
q<-ggplot(data=cvd1, aes(x=date2, y=new_cases_per_million)) +
geom_bar(stat="identity") + labs(x= "Day", y="New Cases, PMP") +
scale_x_date(date_breaks = "1 day") + theme_bw(base_size = 12) +
theme(axis.text.x = element_text(angle = 90))
r<-ggplot(data=cvd1, aes(x=date2, y=total_cases)) +
geom_bar(stat="identity") + labs(x= "Day", y="Cum. Cases, Raw") +
scale_x_date(date_breaks = "1 day") + theme_bw(base_size = 12) +
theme(axis.text.x = element_text(angle = 90))
s<-ggplot(data=cvd1, aes(x=date2, y=total_cases_per_million)) +
geom_bar(stat="identity") + labs(x= "Day", y="Cum. Cases, PMP") +
scale_x_date(date_breaks = "1 day") + theme_bw(base_size = 12) +
theme(axis.text.x = element_text(angle = 90))
ggarrange(p, q, r,s, ncol=2, nrow=2)
dev.off()
str(bl.f)
library(haven)
#install.packages("stats")
require(stats)
#install.packages("graphics")
require(graphics)
#install.packages("ggplot2")
require(ggplot2)
#library("grDevices")
require(grDevices)
#install.packages("gridExtra")
require(gridExtra)
#library("grid")
require(grid)
#install.packages("car")
require(car)
#install.packages("foreign")
require(foreign)
#install.packages("interplot")
require(interplot)
#install.packages("visreg")
#install.packages("devtools")
require(devtools)
#install_github("pbreheny/visreg")
require(visreg)
#install.packages("margins")
require(margins)
#install.packages("coefplot")
#install_github("jaredlander/coefplot")
require(coefplot)
#install.packages("MASS")
require(MASS)
#install.packages("VGAM")
require(VGAM)
#install.packages("stargazer")
require(stargazer)
#install.packages("Hmisc")
require(Hmisc)
#install.packages("Design")
#install.packages("foreach")
#install.packages("psych")
require(psych)
#install.packages("vcd")
require(vcd)
#install.packages("mediation")
#require(mediation)
library(haven)
#library(devtools)
#devtools::install_version("haven", version = "1.1.0")
#install.packages("plyr")
#library(plyr)
library(dplyr)
#install.packages("tidyr")
library(tidyr)
require(multiwayvcov)
require(lmtest)
#devtools::install_github("strengejacke/strengejacke")
#install.packages("sjPlot")
require(sjPlot)
#install.packages("estimatr")
library(estimatr)
#install.packages("fabricatr")
library(fabricatr)
#install.packages("rlang")
require(rlang)
#install.packages("zoo")
library("zoo")
#install.packages("FactoMineR")
library(FactoMineR)
#install.packages("labelled")
library(labelled)
#install.packages("Hmisc")
require(Hmisc)
library(lubridate)
setwd("~/Dropbox/DFID_IOM/Analysis/")
bl_final <- readRDS("bl_final3.rds") # Contacts are updated in bl_final2, participants list amended Feb 15 on bl_final3
cvd.cases <- read.csv("owid-covid-data_gmb_jan2021.csv")
#https://ourworldindata.org/coronavirus/country/gambia
bl.f <- bl_final
bl.f <- copy_labels_from(bl.f,bl_final)
str(bl.f$lost.job)
bl.f <- bl_final
bl.f <- copy_labels_from(bl.f,bl_final)
###################################
#### Employment and Business ######
bl.f$empl_covid <- car::recode(bl.f$empl_covid, "'I did not work, nor did I have my own business.'=0;
'I had a job / worked for someone' = 1;'I was an entrepreneur (self-employed)' = 2")
val_labels(bl.f$empl_covid) <- c('I did not work, nor did I have my own business' = 0, 'I had a job / worked for someone' = 1,
'I was an entrepreneur (self-employed)' = 2)
# Create pre-covid employment dummy: Self employed and employed by someone =1, otherwise 0.
bl.f$wrk.pre.covid <- car::recode(bl.f$empl_covid, "2=1")
bl.f$job <- car::recode(bl.f$job, "2=0")
bl.f$lost.job <- NA
bl.f$lost.job <- ifelse(bl.f$wrk.pre.covid == 0 & bl.f$job == 0, 0,
ifelse(bl.f$wrk.pre.covid == 0 & bl.f$job == 1, 2,
ifelse(bl.f$wrk.pre.covid == 1 & bl.f$job == 0, 1,
ifelse(bl.f$wrk.pre.covid == 1 & bl.f$job == 1, 2,
ifelse(bl.f$job == 1, 2,NA)))))
bl.f$lost.job[which(bl.f$job == 1)] <- 2 # For some reason it's not working the ifelse statement above.
bl.f$lost.job <- labelled(bl.f$lost.job, c('No job before COVID' = 0, 'Lost job after COVID' = 1,
'Kept or got job after COVID' = 2))
str(bl.f$lost.job)
bl.f$lost.job <- as.factor(bl.f$lost.job)
str(bl.f$lost.job)
setwd("~/Dropbox/DFID_IOM/Analysis/II_Covid_Replication")
#save(bl.f, file = "bl.f.RData")
load(file = "bl.f.RData")
bl.f$lost.job <- set_label(bl.f$lost.job, "Lost Job, CV-19. Order: No Job, Lost Job, Kept Job")
var_label(bl.f$lost.job) <-  "Lost Job, CV-19. Order: No Job, Lost Job, Kept Job"
var_label(bl.f$lost.job)
str(bl.f$lost.job1)
str(bl.f$lost.job2)
var_label(bl.f$lost.job) <-  "Lost Job, CV-19. Order: No Job (BL), Lost Job, Kept Job"
var_label(bl.f$lost.job1) <-  "Lost Job, CV-19. Order: Lost Job (BL), No Job, Kept Job"
var_label(bl.f$lost.job2) <-  "Lost Job, CV-19. Order: Kept Job (BL), No Job, Lost Job"
str(bl.f)
str(bl.f$lost.job)
str(bl.f$lost.job1)
str(bl.f$lost.job2)
